function price = binom_option(S0,K,T,r,sigma,N,type)
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%Inputs: S0 -- initial underlying price
%        K -- strike 
%        T -- Maturity
%        r -- interest rate
%        sigme -- volatility
%        N -- number of days
%Output: Erou&Amer&Call&Put
%Method: Binomial Tree
%Example: 
%S0=5;K=10;T=1;r=0.06;sigma=0.3;N=256;
%binom_option(S0,K,T,r,sigma,N)
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
if (nargin < 1)
S0=10;K=10;T=3;r=0.05;sigma=0.2;N=10;type=input('option type :');
%type='AmerCall'or'AmerPut'or'EuroCall'or'EuroPut'
end
dt=T/N; % step size
A=0.5*(exp(-r*dt)+exp((r+sigma^2)*dt));
u=A + sqrt(A^2-1); % appreciation rate
d=1/u;% depreciation rate
q=(exp(r*dt)-d)/(u-d); % risk neutral prob
S = zeros(N+1,N+1); % initialize
euroC = zeros(N+1,N+1); 
amerC = zeros(N+1,N+1);
amerP = zeros(N+1,N+1);
euroP = zeros(N+1,N+1); 


% step I: forward calculation: generate stock path 
for j = 1:N+1
    for i = 1:j
        S(i,j) = S0*u^(j-i)*d^(i-1);
    end
end

% Step II: Backward pricing: calculate call option pricing
switch type
    case 'AmerCall'
        for i = 1:N+1
            amerC(i,N+1) = max(S(i,N+1)-K,0);
        end
        for j = N:-1:1
            for i = 1:j
                amerC(i,j) = max(exp(-r*dt)*(q*amerC(i,j+1)+(1-q)*amerC(i+1,j+1)),S(i,j)-K);
            end
        end
        AmerCall = amerC(1,1)
    case 'AmerPut'
        for i = 1:N+1
            amerP(i,N+1) = max(K-S(i,N+1),0);
        end
        for j = N:-1:1
            for i = 1:j
                amerP(i,j) = max(exp(-r*dt)*(q*amerP(i,j+1)+(1-q)*amerP(i+1,j+1)),K-S(i,j));
            end
        end
        AmerPut = amerP(1,1)
    case 'EuroCall'
        for i = 1:N+1
            euroC(i,N+1) = max(S(i,N+1)-K,0);
        end
        for j = N:-1:1
            for i = 1:j
                euroC(i,j) = exp(-r*dt)*(q*euroC(i,j+1)+(1-q)*euroC(i+1,j+1));
            end
        end
        EuroCall = euroC(1,1)
    case 'EuroPut'
        for i = 1:N+1
            euroP(i,N+1) = max(K-S(i,N+1),0);
        end
        for j = N:-1:1
            for i = 1:j
                euroP(i,j) = exp(-r*dt)*(q*euroP(i,j+1)+(1-q)*euroP(i+1,j+1));
            end
        end
        EuroPut = euroP(1,1)
end

end